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NA44 uses a 512 channel Si pad array covering 1.5 < rj < 3.3 to study charged hadron production 
in 158 A GeV Pb + Pb collisions at the CERN SPS. We apply a multiresolution analysis, based on 
a Discrete Wavelet Transformation, to probe the texture of particle distributions event-by-event, 
allowing simultaneous localization of features in space and scale. Scanning a broad range of mul- 
tiplicities, we search for signals of clustering and of critical behavior in the power spectra of local 
density fluctuations. The data are compared with detailed simulations of detector response, using 
heavy ion event generators, and with a reference sample created via event mixing. An upper limit 
is set on the probability and magnitude of dynamical fluctuations. 
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The main experimental challenge in relativistic heavy 
ion collisions is to find evidence for the expected QCD 
phase transition at high temperature. Deconfinement 
and chiral symmetry restoration are expected to take 
place during the hot, strongly interacting stage early in 
the collision. As a phase transition in such collisions is 
inherently a multiparticle phenomenon, multiparticle ob- 
servables, defined on event-by-event basis, are of great 
interest. Recently published event-by-event analyses of 
the 158 GeV/ A Pb+Pb data either analyze a small num- 
ber of events H in great detail, or analyze properties of a 
large ensemble of events using a single observable (jpt) to 
compare different ensemble averages pi. In the first case, 
accumulation of feature information from large data sets 
remains an open issue. In the second case, an ensemble 
average on a set of post- freeze- out events is not represen- 



tative of the pre- freeze- out history of those events, due to 
the dramatic non-stationarity of the open system, with 
a consequent lack of ergodicity. Violations of ergodicity 
generally happen in the course of phase transitions |3j. 

We concentrate on texture, or local fluctuation observ- 
ables, analyzing single events independently to determine 
the scale composition of fluctuations. In the following, 
we may omit the term "local" , but we will always talk 
about fluctuations in the particle density from one point 
to another within a single event, i.e. in the local sense, 
as opposed to fluctuations of global quantities from one 
event to another. 

In 1985, L. Van Hove formulated a model of quark- 
gluon plasma hadronizationH with a first order phase 
transition. Longitudinal expansion of the colliding sys- 
tem, with particle formation via string or color flux tube 
breaking, can result in plasma droplets as large as a few 
fm across. The droplets hadronize by deflagration |5j] . 
This is expected to result in dN / dy distributions with 
bumps or spikes on top of an otherwise smooth struc- 
ture. Other models H also predict bubbles of one phase 
embedded in the other. 

In the absence of a direct, event-by-event observable- 
based test of these predictions, the picture had been fur- 
ther developed |7], || in order to connect it with the tra- 
ditional observables such as the Tot slope parameter T 
and the baryon and strangeness chemical potentials: the 
hadron "temperatures" T in the SPS data are higher than 
lattice QCD predictions for a phase transition tempera- 
ture. Using a first order phase transition hydrodynamical 
model with a sharp front between the phases, Bilic et al. 
[M M concluded that a QGP supercooling and hadron gas 
superheating is a consequence of the continuity equations 
and of the requirement that the entropy be increased in 
the transition. In the case of bubbles in the QGP phase, 
the plasma deflagrates; otherwise, it detonates. A direct 
measurement of the hadron texture at freeze-out, if it de- 
tects presence of the droplets/bubbles, could provide an 
argument in favor of the first order phase transition. 

The order of the confinement phase transition is still 
under debate. It is a fluctuation driven first order tran- 
sition [J9], [lCJ] in SU(3) with three massless quarks, but 
second order in the case of finite mass |ll| or infinitely 
massive || [Hj strange quarks. A tricitical point may 
exist, separating the first order transition from a second 
order transition with the same critical exponents as the 
3D Ising model ||. For a second order phase transi- 
tion, local fluctuations of isospin or enhanced correlation 
lengths may be observable |1J, O. Large scale corre- 
lations formed early in the collision are more likely to 
survive diffusion in the later stages. Small scale fluctu- 
ations, on the contrary, are more easily washed out by 
diffusion due to secondary scattering among the hadrons 
JL4|. Consequently an analysis method which can iden- 
tify fluctuations on any scale is desirable. In this paper, 
we utilize a Discrete Wavelet Transformation, which has 
this property. 

A wavelet is a function, everywhere except for a well 



localized spot. For a pad detector, the discrete positions 
of the spots correspond naturally to the pad positions, 
and the possible scales are multiples of the pad sizes. The 
scale is an analog of a Fourier frequency. Location has no 
analog in the Fourier transform, and it provides an addi- 
tional degree of analytical power, which explains much of 
the success that wavelets met in the field of data process- 
ing and pattern recognition. (Examples of Fourier-based 
analyses of large scale azimuthal texture in the field of 
relativistic heavy ion collisions exist as well (lj, [f(| ; this 
is how the elliptic flow at relativistic energies was mea- 
sured.) The binning of charged particle density inherent 
in measurements with a segmented detector such as a Si 
pad detector makes the Haar wavelet a natural choice of 
analyzing function; a Haar wavelet is a step function with 
given width, oscillating around zero with a single period. 
Discrete Wavelet Transformation (DWT)|l7| quanti- 
fies contributions of different <f> and r\ scales to the event 
texture. We use DWT to test for possible large scale en- 
hancement, as a function of the collision centrality. We 
report the DWT power spectrum in pseudorapidity r\ and 
azimuthal angle 0, for different charged particle multi- 
plicities. We use mixed events to remove trivial fluctua- 
tions and background effects. 



II. EXPERIMENTAL SETUP 

The experimental setup |18[ is shown in FigJjL The 
NA44 Si pad array installed 10 cm downstream from the 
target, in the magnetic field of the first dipole, measured 
ionization energy loss of charged particles in its 512 300 
/xm thick Si pads. The plastic scintillator TO (two rect- 
angles seen in Fig.nT) was used for a centrality trigger. 
The SPS beam was collimated to a 1 x 2 mm profile. 
TO covered 1.4 < i] < 3.7 for an 77-dependent fraction 
of azimuthal angle, 0.22 < A<j>/2tt < 0.84 respectively. 
The silicon detector had inner radius 7.2 mm and outer 
radius 43 mm, covering 1.5 < r\ < 3.3. The detector 
was split radially into 16 rings of equal rj coverage. Each 
ring was further divided azimuthally into 32 sectors of 
equal angular coverage to form pads. The pads were 
read out by AMPLEX |19[ chips, one chip per sector. 
(5-electrons, produced by the Pb beam traversing the tar- 
get, were swept away to one side by the dipole magnetic 
field (< 1.6 T). Only the <5-electron-free side was used in 
this analysis. Only 4 of the remaining 256 channels were 
inoperative. 

An amplitude distribution from a typical channel, ob- 
served in the physics run and digitized with a 256 channel 
ADC is shown on FigM Channel pedestals had, on the 
average, FWHM = 0.48 < dE > of 1 MIP. In the tex- 
ture analysis, every event was represented by a 2D array 
of the calibrated digitized amplitudes of the channels ( 
an amplitude array). Empty target runs were used to 
measure the background, and cross-talk in the detector 
was evaluated off-line. 



a) 





FIG. 1: a) The experimental setup: the target, the Si pad 
array and the TO scintillation counter. See text for a descrip- 
tion of the detectors, b) The setup exposed to an RQMD 
event (GEANT simulation). Magnetic field is on. 
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FIG. 2: Digitized amplitude distribution from channel 1 of the 
Si pad array. The smooth curve shows a minimum x 2 Landau 
fit performed in course of the amplitude calibration. The 
pedestal, the single and double hit peaks are distinguishable. 



<fi = 4>{rf , (/)'). This makes the observable multiplicity 
distribution d 2 N / d<f>' dr\' (in the presumed coordinates) 
differ from a simple function of r/': 



d 2 N 1 dN 

d<t>' drf ^ 2^"*/ 



(1) 



III. ANALYSIS TECHNIQUE 
A. Detector calibration 

The NA44 spectrometer information was not used in 
this analysis, which focussed on the Si pad array data. 

ADC pedestals were fitted channel by channel with a 
realistic functional shape, determined from low multiplic- 
ity events in a minimum bias triggered run. Amplitude 
calibration of the Si detector was carried out channel by 
channel, by fitting the amplitude distribution with a sum 
of single, double, triple, etc. (up to septuple) minimum 
ionizing particle Landau distributions [2TJ with variable 
weights. The Landau distributions were numerically con- 
voluted with the pedestal shape to account for noise in 
the fit. A typical fit from a single channel is shown in 
Fig. g. Parameters of the fit were used to simulate noise 
in a GEANT-based detector response Monte Carlo code. 

An offset of the event vertex with respect to the de- 
tector's symmetry axis results in a non-trivial functional 
dependence between the actual r\ and 0, and the r/, <f>' 
presumed based on the "ideal" geometry: rj = rj(r]' , </>'), 



In the true coordinates rj and 4>, the inequality [D becomes 
an equality. However, the detector's acceptance area in 
the true coordinates becomes distorted. In the following 
we will refer to this as a "Jacobian effect" . The Jacobian 
effect, obviously, contributes to the event textures, espe- 
cially on the large scale, and needs to be evaluated and 
corrected for. 

From Eq.pl, the criterion of the true coordinate basis 
(rj, (f>) emerges naturally: it is the basis which makes the 
observable d 2 N/d(j)drj independent of <i>. The minimiza- 



tion problem was solved numerically with MINUIT|21|, 
and the resulting offsets are within the tolerance of the 
detector/beam position. Cross-talk between the electron- 
ics channels is a detector-related correlation phenomenon 
and introduces a "texture" effect of its own. Both global 
cross-talk in the AMPLEX read-out chip |]19[ and read- 
out board cross-talk are expected. In our detector with 
512 channels, there are 512 x (512 - l)/2 = 130816 two 
channel pairs (unordered), all of which were subjected to 
covariance analysis off-line. To magnify the non-trivial 
instrumental contribution to the covariance matrix ele- 
ments, we analyzed covariances not between the ampli- 



tudes Aa of channels i themselves, but between 
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Otherwise, the dominant contributor to the cov(Ai,Aj) 
is the trivial variation of the event's common multiplic- 
ity [p6| . Using this method, we concluded that the effec- 
tive cross-talk coupling was non-negligible only for neigh- 
boring channels within the same chip; it was found to 
be 8.5%. As a remedy, a chip-wise (i.e. sector-wise) 
event mixing technique including cross-talk in the ref- 
erence sample was used to construct a reference event 
sample. 
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FIG. 3: Double differential multiplicity distributions of 
charged particles plotted as a function of azimuthal angle </> 
(with different symbols representing different rings) and of 
pseudorapidity 77 (with different symbols representing differ- 
ent sectors). The <j> and 77 are in the aligned coordinates. 

The double differential multiplicity data (Fig. H) illus- 
trate the quality of the detector operation, calibrations, 
geometrical alignment and Jacobian correction. The data 
set is composed of two pieces, obtained by switching the 
magnetic held polarity: a negative polarity run is used 
for sectors 9 to 24 (range of n/2 < 4> < 3ir/2); a posi- 
tive polarity run is used for sectors 1 to 8 and 25 to 32 
(range of < <fi < tt/2 and 37r/2 < <j> < 2ir). The reason 
to disregard one side of the detector is additional oc- 
cupancy due to (5-electrons, as was explained in section 
II. Figure |3j demonstrates the quality of alignment as 
well, since the rj and <f> along the horizontal axes are the 
aligned coordinates. Any geometrical offset of the detec- 
tor makes acceptances of different pads non-equal and de- 
pendent on the pad position. The acceptance of each pad 
has been calculated in the aligned coordinates, and the 
d 2 N J d(j) dij uses the actual acceptances dcj). The shape 
of the 4> dependence of d 2 N/ d(j)drj (left panel of Fig. ||) 
is flat as it should be for an event ensemble with no re- 
action plane selection. The 77 dependence (right panel 



of Fig. 0) shows increasing multiplicity towards midra- 
pidity, as is expected. As can be seen from Fig.g, the 
detector's acceptance is asymmetric around midrapidity. 
A correction for the cross-talk has been applied. 



B. Discrete Wavelet Transformation (DWT) 

Discrete wavelets are a set of functions, each having 
a proper width, or scale, and a proper location so that 
the function differs from only within that width and 
around that location. The set of possible scales and lo- 
cations is discrete. The DWT formalizes the two dimen- 
sional particle distribution in each Pb + Pb collision in 
pseudorapidity 77 and azimuthal angle (j) by performing 
an image analysis - transforming the event into a set of 
functions orthogonal with respect to scale and location 
in the (77, </>) space. We accumulate texture information 
by averaging the power spectra of many events. 

The simplest DWT basis is the Haar wavelet, built 
upon the scaling functionp^ g(x) — 1 for < x < 1 and 
otherwise. The function 
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is the wavelet function |38|. 

If the interaction vertex lies on the detector's sym- 
metry axis, every pad's acceptance is a rectangle in 
the (0, 77) space. Then, the Haar basis is the natural 
choice, as its scaling function in two dimensions (2D) 
G(4>, 77) = g(<fi)g(r]) is just a pad's acceptance (modulo 
units). We set up a two dimensional (2D) wavelet basis: 
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The scaling function in 2D is G((j>, rj) — g{<p)g(rj). As 
in Eqjj, we construct G m: ij((j),r]) where m is the integer 
scale fineness index, and i and j index the positions of bin 
centers in <f> and 77 (1 < m < 4 and 1 < i, j < 16 because 
we use 16 = 2 4 rings and 16 sectors). Different values 
of A (denoted as 0, 77, and tprj) distinguish, respectively, 
functions with azimuthal, pseudorapidity, and diagonal 
texture sensitivity: 



p4 



f(4)g(v), F n = gMttv), F* 1 = f{4>)Kri) (5) 



Then, F^ i ■ with integer m, i, and j are known ||17| to 
form an orthonormal basis in the space of all measur- 
able functions defined on the continuum of real numbers 
L 2 (R) . Fig. ffl shows the wavelet basis functions F in two 
dimensions. At first glance it might seem surprising that, 
unlike the ID case, both / and g enter the wavelet basis 
in 2D. Fig. [| clarifies this: in order to fully encode an 
arbitrary shape of a measurable 2D function, one con- 
siders it as an addition of a change along <j> (f(4>)g(rf), 
panel (b)), a change along 77 (g(4>)f{r]), panel (c)), and 




FIG. 4: Haar wavelet basis in two dimensions. The three 
modes of directional sensitivity are: a) diagonal b) azimuthal 
c) pseudorapidity. For the finest scale used, the acceptance of 
a Si pad would correspond to the white rectangle drawn "on 
top" of the function in panel a). Every subsequent coarser 
scale is obtained by expanding the functions of the previous 
scale by a factor of 2 in both dimensions. 



a saddle-point pattern (f(<f>)f(r]), panel (c)), added with 
appropriate weight (positive, negative or zero), for a va- 
riety of scales. The finest scale available is determined 
by the detector segmentation, while the coarser scales 
correspond to successively rebinning the track distribu- 
tion. The analysis is best visualized by considering the 
scaling function G mj ij((f>,rj) as binning the track distri- 
bution p(<fi,rj) in bins i,j of fineness to, while the set of 
wavelet functions F^ i t j(<f>, rj) (or, to be exact, the wavelet 
expansion coefficients (p, F^ i ■)) gives the difference dis- 
tribution between the data binned with given coarseness 
and that with binning one step finer. 

While the DWT analyzes the object (an image, a se- 
quence of data points, a data array) by transforming it, 
the full information content inherent in the object is pre- 
served in the transformation. 

We adopt the existing [£2| ID DWT power spectrum 
analysis technique and expand it to 2D. The track density 
in an individual event is p{4>, rf) and its local fluctuation 
in a given event is a 1 = (p — p, p — p) , where p is the 
average p (over the acceptance) in the given event |3S[ . 

Using completeness of the basis, we expand 
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Notice that p, being constant within the detector's 
rectangular acceptance, is orthogonal to any -F^* j 
with to > I. Due to the orthonormality condition 
( F m,i,j> F m>,i',j>) = d\,\'Sm,m'5i,i'Sj,j', the p - p compo- 
nents for different scales do not form cross-terms in the a 2 
sum, and the sum contains no cross-terms between p and 
p for the four observable scales. Instead of a (p, G m= 5^j) 
set, the Si detector energy amplitude array - its closest 
experimentally achievable approximation - is used as the 
DWT input. We used the WAILI [§| software library to 
obtain the wavelet decompositions. 

The Fourier power spectrum of a random white noise 
field is known to be independent of frequency [|i| . We 
are looking for dynamical textures in the data, and there- 
fore would like to treat the random white noise case as 
a "trivial" one to compare with. Therefore it is interest- 
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FIG. 5: Understanding the analyzing potency of the DWT 
power spectra: a) for a checkerboard pattern b) for a smooth 
gradient pattern c) for a sample of a thousand random white 
noise images - in this case the average power spectrum is 
shown. 



ing to reformulate this property for wavelets, where scale 
plays the same role as frequency in Fourier analysis. 

To do that, we link scales with frequencies, or in other 
words, we must understand the frequency spectra of the 
wavelets. The Fourier images of ID wavelet functions oc- 
cupy a set of wave numbers whose characteristic broad- 
ness grows with scale fineness m as 2 m ; 2 2m should be 
used in the 2D case. Discrete wavelets of higher or- 
ders have better frequency localization than the Haar 
wavelets. Despite this advantage, we use Haar because 
only Haar allows one to say that the act of data taking 
with the (binned !) detector constitutes the first stage of 
the wavelet transformation. 

In 2D, we find it most informative to present the three 
modes of a power spectrum with different directions of 
sensitivity P^ v (m), P^(m), P r '(m) separately. We define 
the power spectrum as 



v H = iE<ftC> 2 - 



(7) 



where the denominator gives the meaning of spectral den- 
sity to the observable. So defined, the P x (m) of a random 
white noise field is independent of to. 

In order to illustrate the sensitivity of the wavelet 
transformation to texture features of the different scales, 
we have applied the wavelet transform to three test pat- 
terns, shown in Fig. 0. All patterns are 16 x 16 pixel 
matrices. The left hand side shows the test pattern, and 
the right shows the power spectrum resulting from the 
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FIG. 6: Power spectra of 7 x 10 3 events in the multiplicity 
bin 326 < dN/ dn < 398 (between w 6% and 10% centrality). 
O _ true events, A - mixed events, □ - the average event. 



wavelet transform. Pattern a), a checkerboard, has struc- 
ture only on the finest scale and all power components 
of scales coarser than 4 are zero. Pattern b) has exactly 
the opposite scale composition; the slow gradation be- 
tween black and white corresponds to a structure on the 
coarsest scale, as seen in the accompanying power spec- 
trum. Smoothness of the gradient means that neighbor- 
to-neighbor changes do not add much to the pattern once 
the overall trend (the large scale feature) is taken into ac- 
count. 

These two examples illustrate the property of scale lo- 
calization, made possible by virtue of the scale orthog- 
onality of the basis. Patterns encountered in multiple 
hadron production involve a variety of scales, and yet 
they are more likely to be of type b), rather than a). An 
important conclusion follows immediately: in this type 
of measurement, large acceptance, like the one used in 
this analysis, rather than fine segmentation, is the way 
to accomplish sensitivity. 

Case c) shows patterns that arise from white noise. 
They produce signals in the power spectrum indepen- 
dent of scale, as expected. In the first approximation, the 
white noise example provides a base-line case for compar- 
isons in a search for non-trivial effects. 

Figure |6| shows the power spectra measured in Pb + Pb 
for one multiplicity range. The unit on the vertical scale 
(a 2 /{dEMip) 2 ) is chosen so that the power of fluctua- 
tions whose variance a 1 equals the squared mean energy 
loss by a minimum ionizing particle traversing the de- 
tector, is the unit. The first striking feature is that the 
power spectra of physical events are indeed enhanced on 
the coarse scale. The task of the analysis is to quan- 
tify and, as much as possible, eliminate "trivial" and 
experiment-specific reasons for this enhancement. 



duced texture as statistical fluctuations cancel (shown as 
□ in FigJI) . Average events retain the texture associated 
with the shape of d 2 N / dip dr), with the dead channels and 
the finite beam geometrical cross-section (though this is 
only partially visible in the average event, due to the fact 
that event averaging is done without attempting to se- 
lect events according to the vertex position). P x (m) is 
proportional to the variance, or squared fluctuation a 2 . 
Therefore, for Poissonian statistics of hits in a pad, the 
event averaging over M events should decrease P x (m) by 
a factor of M. The average event whose power spectrum 
is shown on Fig. g] is formed by adding 7 x 10 3 events, 
however its P x (m) is down less than 7 x 10 3 compared 
to that of the single events. This demonstrates that the 
average event's texture is not due to statistical fluctua- 
tions, but rather, predominantly due to the systematic 
uncertainties listed. Consequently, we can use the aver- 
age event's P x (m) to estimate the magnitude of the static 
texture-related systematics. As seen from Fig. g, the sys- 
tematics are far below the P x (m) of single events (true 
or mixed), with the exception of pseudorapidity, where 
non-constancy of dN/ dij over the detector's acceptance 
is visible. 

The way to get rid of the "trivial" or static texture 
is to use mixed events, taking different channels from 
different events. The mixed events preserve the texture 
associated with the detector position offset, the inherent 
dN/ dn shape and the dead channels. This is static tex- 
ture as it produces the same pattern event after event 
while we are searching for evidence of dynamic texture. 
We reduce sources of the static texture in the power spec- 
tra by empty target subtraction and by subtraction of 
mixed events power spectra, thus obtaining the dynamic 
texture P x (m) true — P x (m) m j x . In order to reproduce the 
electronics cross-talk effects in the mixed event sample, 
the mixing is done sector-wise, i.e. the sectors constitute 
the subevents subjected to the event number scrambling. 

We continue with a brief summary of the systematic 
errors in the measurements of the DWT dynamic tex- 
ture observable Ptrue — Pmix- Static texture and dy- 
namic background texture present the largest problem in 
the search for the phase transition-related dynamic tex- 
ture via power spectra of local fluctuations. The method 
of solving the problem is comparison with the reference 
sample created by event mixing. Thus the Ptrue — Pmix 
observable was created. For comparison with models, a 
Monte Carlo simulation of the Si detector is used. It in- 
cludes the known static texture effects and undergoes the 
same procedure to remove the effects. The "irreducible 
remainder" is the residual effect which may 

1. survive the elimination procedure 



C. Identification and control of systematic errors 

The average event, formed by summing amplitude im- 
ages of the measured events in a given multiplicity bin, 
and dividing by the number of events, has a much re- 



2. emerge as a difference between the data, subjected 
to the elimination procedure, and the MC analyzed 
in the same manner. 

Table m lists the sources of static texture and summarizes 
the methods of their treatment. We group the background 



texture sources according to similarity of manifestation 
and treatment, into 

• statistical fluctuations 

• static texture 

• background dynamic texture 

The statistical fluctuation is the most trivial item in this 
list. Both event mixing (provided that mixing is done 
within the proper multiplicity class) and MC comparison 
solve this problem. The statistical fluctuations do not 
result in irreducible systematic errors. 
The static texture group includes: 

• geometrical offset of the detector with respect to 
the beam's "center of gravity" in the vertical plane 

• dead pads 

• dN / dr\ shape - a genuine large scale multiparticle 
correlation sensitive to the physics of the early stage 
of the collision 

Cleanliness of the static texture elimination via event 
mixing has been checked by simulating the contribut- 
ing effects separately. First, by running the detector re- 
sponse simulation on MC-generated events without the 
beam/detector offset and with a beam of thickness it 
was ascertained that the remaining dynamic texture is 
very small compared with the systematic errors due to 
the background Si hits and the beam geometrical cross- 
section, for all scales and all directional modes A. Due 
to the finite size of the multiplicity bin, the mixed events 
consist of subevents coming from events of different total 
multiplicity. With the sector-wise mixing, this causes an 
additional sector-to-sector variation of amplitude in the 
mixed events, thus resulting in an enhancement of P„ lix 
primarily on the finest scale, with respect to Pf rue - On 
Fig. |t], this effect can be seen as the P{ rue — P mix values 
progressively grow negative with multiplicity in the finest 
scale plot. However, as can be seen on the same figure, 
the effect is small compared with the total systematic 
error bars shown as boxes. 

The background dynamic texture group includes: 

• elliptic and directed flow 

• finiteness of the beam cross-section 

• background hits in the Si 

• channel-to-channel cross-talk 

Elliptic and directed flow, observed at SPS [16[, are 
large scale dynamic texture phenomena of primarily az- 
imuthal (elliptic) and diagonal (directed flow) modes. 
Because both reaction plane and direction angle vary 
event by event, the respective dynamic textures can not 
be subtracted by event mixing, unless the events are clas- 
sified according to their reaction plane orientation and 



the direction angle, with mixing and Ptrue — Pmix sub- 
traction done within those classes. Neither reaction plane 
nor direction angle was reconstructed in the present anal- 
ysis, and the Ptrue — Pmix (especially that of the az- 
imuthal and diagonal modes on the coarse scale) retain 
the elliptic/directed flow contribution. The effects of flow 
on dynamic texture observables are smaller than other 
texture effects, so they can not be singled out and quan- 
tified in this analysis. 

The finite beam cross-section effect belongs to this 
group, despite the fact that a very similar effect of geo- 
metrical detector/beam offset has been classified as static 
texture. An effect must survive mixing with its strength 
unaltered in order to be fully subtracted via event mix- 
ing. Preserving the effect of the random variations in the 
Pb + Pb vertex on the power spectra in the mixed events 
requires classification of events according to the vertex 
position and mixing only within such classes. This re- 
quires knowledge of the vertex for each event, which is 
not available in this experiment. Therefore, MC simula- 
tion of the beam profile remains the only way to quantify 
false texture arising from vertex variations. MC studies 
with event generators show that the beam spatial ex- 
tent and the resulting vertex variation is the source of 
the growth of the coarse scale azimuthal texture corre- 
lation with multiplicity (see Fig. M). Uncertainty in our 
knowledge of the beam's geometrical cross-section must 
be propagated into a systematic error on Ptrue — Pmix- 

The other two effects in this group are difficult to sep- 
arate and simulate and the error estimate reflects the 
combined effect. The systematic errors were evaluated 
by removing the Pb target and switching magnetic field 
polarity to expose the given side of the detector to 8- 
electrons (from the air and TO), while minimizing nu- 
clear interactions. This gives an "analog" generator of 
uncorrelated noise. All correlations (i.e. deviations of 
P X (iTi)true from P x (rn) m i X ) in this noise generator are 
treated as systematic uncertainties. Thus this compo- 
nent of the systematic error gets a sign, and the sys- 
tematic errors are asymmetric. The effect of increasing 
texture correlation (for diagonal and azimuthal modes) 
with multiplicity on the coarse scale, attributed to the ge- 
ometrical offset of the detector with respect to the beam 
(the leading one in the static group), is present in the 
switched polarity empty target runs as well. For this 
reason, it was impossible to disentangle the background 
dynamic contribution on the coarsest scale. In Table |, 
the "irreducible remainder estimate" for the diagonal, 
coarse scale is bracketed with two numbers, which form 
the lower and upper estimates. The lower estimate is ob- 
tained by taking the scale one unit finer and quoting its 
number. This, indeed, sets the lower limit because the 
deviations of P x (m) t rue from P x (m) m i X generally grow 
with scale coarseness. The upper limit is set by ascribing 
the entire texture correlation, observed in the <5-electron 
data, to the background hits and channel cross-talk, and 
ignoring the fact that significant portion of it must be 
due to the vertex fluctuation (finite beam profile). This 
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FIG. 7: Multiplicity dependence of the texture correlation. 
O - the NA44 data, • - RQMD. The boxes show the system- 
atic errors vertically and the boundaries of the multiplicity 
bins horizontally; the statistical errors are indicated by the 
vertical bars on the points. The rows correspond to the scale 
fineness m, the columns - to the directional mode A (which 
can be diagonal (pr/, azimuthal <f>, and pseudorapidity rfj. 



upper limit is likely to be a gross overestimation, and in 
Fig. U\ we show systematic errors, obtained by adding in 
quadrature the finite beam error with the background hit 
error. 



IV. RESULTS 

Fig. presents a comparison of the DWT dynamic tex- 
ture in the measured and RQMD-simulated pq| Pb + Pb 
collision events. The three directional sensitivity modes 
(diagonal 4>rj, azimuthal 0, and pseudorapidity rj) have 
four scales each, so that there are 12 sets of points in 
the DWT dynamic texture as a function of the charged 
multiplicity dN c h/ drj bin. The systematic errors on the 
points (shown by vertical bars) have been evaluated fol- 
lowing the procedure described in detail in Section [II. 

Fig.y demonstrated that the major fraction of the ob- 
served texture exists also in mixed events. A detailed 
account of the causes was discussed in the preceding sec- 
tion, including known physics as well as instrumental ef- 
fects. It is therefore clear that the observable most di- 
rectly related to the dynamical correlations/fluctuations 
is not P x (m), but P x (m)true ~ P X ('m)mix- This quan- 
tity, normalized to the RMS fluctuation of P x (m) m i X , 
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FIG. 8: Confidence coefficient as a function of the fluctuation 
strength. RMS mtx denotes ^/{P x {l) 2 mix - (P A (1)™*} 2 ). 
The multiplicity bin is 326 < dN/ dn < 398 (6-10% central- 
ity), as in Fig.f 



can be used to characterize the relative strength of local 
fluctuations in an event. The distribution for different A 
(or directions) is plotted on Figure H in an integral way, 
i.e. as an a(x) graph where for every x, a is the fraction 
of the distribution above x. 
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where £ denotes the fluctuation strength 
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RMS{P X (l\ 



(8) 
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and dN/ dl; is the statistical distribution of £, obtained 
from the experimentally known distributions of P x (l)true 
and P x (l) m ix- Expression p] is constructed to be sensi- 
tive to the difference between P x (l)true and P x (l) m ix, 
while minimizing detector specifics to enable comparison 
between different experiments in future. The latter is 
accomplished by normalizing to RMS m i X - This normal- 
ization also eliminates the trivial multiplicity dependence 
of the observable. 

The fluctuation strength observable provides a limit 
on the frequency and strength of the fluctuations and 
expresses the result in a model-independent way. The 
confidence level with which local fluctuations of a given 
strength (expressed through the event by event observ- 
ables via Eq. ^) can be excluded is then 1 — a. Fluc- 
tuations greater than 3 x RMS m i x are excluded in the 
azimuthal and pseudorapidity modes with 90% and 95% 
confidence, respectively. The monotonic fall of the curve 
is consistent with the absence of abnormal subsamples in 
the data. 

RQMD events were fed into the GEANT detector re- 
sponse simulation and analyzed using the same off-line 
procedure as used for the experimental data. The de- 
tector offset with respect to the beam center of gravity 
and the beam profile were included in the simulation. In 
a separate simulation run, the beam profile was identi- 
fied as the cause of the rise of the azimuthal dynamic 



texture with the multiplicity on the coarse scale. In our 
experiment, this purely instrumental effect dominates the 
azimuthal component of the DWT dynamic texture. 

The most apparent conclusion from Fig. is that a 
large fraction of the texture (seen on Fig. ||) is not dy- 
namic i.e. not different between true and mixed events. 
Being monotonic (or absent), the change of the data 
points with multiplicity does not reveal any evidence of 
a region of impact parameters/baryochemical potentials 
with qualitatively different properties, such as those of 
a critical point neighborhood. The RQMD comparison 
confirms that particle production via hadronic multiple 
scattering, following string decays (without critical phe- 
nomena or phase transition) can explain the observed re- 
sults when detector imperfections are taken into account. 
More detailed discussion of the implications of these data 
on various phase transition models will be given in Sec- 



tion VI 



V. SENSITIVITY 

Interesting physics can manifest itself in the ensemble 
probability density distributions as well as in the event- 
by-event (EbyE for short) observables. To illustrate the 
power of the EbyE observable we used, we should con- 
struct final states of charged particles indistinguishable 
from the point of view of "traditional" , or ensemble- wise 
observables, such as 

1. dN/ dy distribution 

2. dN J dpT, 1/mTdN/drriT distribution etc. 

3. multiplicity distribution 

and compare the sensitivity of the above-mentioned ob- 
servables with that of the EbyE one. 

A sensitivity study was performed using a multifireball 
event generator created specially for this purpose. The 
generator produces textures of known magnitude by sim- 
ulating the observed multiplicity as arising from an ar- 
bitrary number of fireballs. Correlations among groups 
of particles arise when the particles come from the same 
fireball. We do not suggest that the physics of Pb + Pb 
collisions is properly described by a superposition of fire- 
balls of a fixed size. Rather, we use the fireballs as a way 
to generate controlled multiparticle correlations. 

This picture is inspired by Van Hove's scenario Q 
of a first order phase transition via droplet fragmenta- 
tion of a QGP fluid. We measure texture in two direc- 
tions, spanned by polar and azimuthal angles, and are 
also sensitive to the spatial texture of longitudinal flow. 
For boost-invariant expansion p6| two droplets, separated 
along the longitudinal coordinate, will be separated in y 
and X]. As long as there is longitudinal expansion, a spa- 
tial texture will be manifested as (pseudo)rapidity tex- 
ture. In the multifireball event generator, we generate 
the pseudorapidity texture explicitly, omitting the spa- 
tial formulation of the problem. The total pt of each 



400 
200 



"S 400 
D) 
a- 200 

P 200 



j^L 



* -A-A- 



O^OoO 1 



i0 o 



*■> 



I , .D^ 



^o. .o-%. 



>A-A 



00 



OQ 



^ncw 



^ 



00 



o o 



?°oo^ I 



300 

200 

Z 100 

"° 

300 

200 

100 



12 3 4 5 6 

y 

FIG. 9: dN/ dy distribution of charged particles in the multi- 
fireball event generator in four individual events with different 
number of fireballs: A - 2 fireballs, □ — 4 fireballs, O - 8 fire- 
balls, O ~ 1 6 fireballs. One can see how the texture becomes 
smoother as the number of fireballs increases. We remind the 
reader that the detector's active area covers 7r azimuthally and 
pseudorapidity 1.5 to 3.3. In general, acceptance limitations 
make it more difficult to detect dynamic textures. 



fireball is 0; its total pz is chosen to reproduce the ob- 
served dN/ dy of charged particles by Lorentz-boosting 
the fireballs along the Z direction, keeping the total p 
of an event at in the rest frame of the colliding nuclei. 
The fireballs hadronize independently into charged and 
neutral pions and kaons mixed in a realistic proportion. 
By varying number of particles N p per fireball, one varies 
"grain coarseness" of the event texture in r\. 

To illustrate the discussion, Fig.pl presents examples of 
dN / dy distributions in four events with different number 
of fireballs. The dynamic textures seen on the figures 
are peculiar to these particular events and are gone after 
dN / dy of many events are added. 

We simulated average fireball multiplicities of 10, 50, 
90 (with RMS fluctuation of 3) and larger. Fig. [l(] shows 
comparison of our data with the simulated pseudorapid- 
ity texture. With ~ 10 4 events, the detector+software 
can differentiate between the cases of 50 and 90 parti- 
cle fireballs. The signal grows with the charged particle 
multiplicity and with N p . Fig. |l0| provides quantitative 
information on the sensitivity of the texture measure- 
ments by relating the expected strength of response to 
the strength of texture via Monte Carlo simulation. The 
sensitivity is limited by systematic errors of the mea- 
surement, discussed in Section 111. Nevertheless, it is 



instructive to compare sensitivity of this method with 
other methods; in particular with two point correlators. 
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The sensitivity of the method is remarkable indeed if 
one takes into account that statistics in the fifth multi- 
plicity bin for each of the three event generator points is 
below 3 x 10 4 events - too scarce, e.g., to extract three 
source radius parameters via HBT analysis even with a 
well optimized spectrometer! 

The use of two particle correlation in rapidity i?2(y) 
to search for droplets was discussed for pp collisions at 
s/s = 1.8 TcV (at FNAL)J||l. R 2 was reported to de- 
crease with multiplicity, so that it would not be expected 
to be visible for dN / dy above «20; the signal would be 
weaker in a scenario with correlated droplets. In con- 
trast, the wavelet transformation retains sensitivity at 
high multiplicity, as we see in FigJI^. In the fifth mul- 
tiplicity bin, with total number of hadrons at freeze-out 
around 1.5 x 10 2 , a typical fraction of particles coming 
from the same fireball for the clustering parameters of 
50 (90) would be 3% (6%) go). In either case there is 
little hope of seeing any trace of such dynamics either in 
ensemble-averaged dN / dy or in dN / dy of a single event, 
but the systematic difference between the power spectra 
of the real and mixed events, integrated over multiple 
events, nevertheless reveals the difference. The data are 
consistent with clustering among < 3% of the particles. 



VI. DISCUSSION 

The order of the expected QCD phase transitions is 
known to be a complex issue for realistic current masses 
of quarks in the system of a finite size. It is generally ex- 



pected that a first order phase transition would be easier 
to observe. Our dynamic texture measurement tests the 
hypothesis of the first order phase transition via QGP 
droplet hadronization^ in a way more direct than in- 
terpretation of pt spectra involving latent heat. Our re- 
sult can be used to constrain phenomenological quantities 
which represent basic QCD properties and affect texture 
formation in this class of hadronization models ||, g, |) . 
Such quantities are the energy flux, or rate at which the 
QGP transmits its energy to hadrons [£9[ |3(| , critical size 
of the QGP droplet Q, and initial upper energy density 
of the transition e' . 

The specific experimental signature of second order 
phase transition (known since the discovery of critical 
opalescence [glj) is the emergence of critical fluctuations 
of the order parameter with an enormous increase of the 
correlation lengths. However, for physical quark masses 
Rajagopal and Wilczek |lj, |32[ argued that due to close- 
ness of the pion mass to the critical temperature, it 
would be unlikely for the correlation volumes to include 
large numbers of pions, if the cooling of the plasma and 
hadronization proceeds in an equilibrated manner. If, 
on the contrary, the high temperature configuration sud- 
denly finds itself at a low temperature, a self-organized 
criticality regime settles in, and the critical local fluctu- 
ations develop fully (T|, ||. 

The NA44 data reported here signifies absence of dy- 
namical fluctuations on the scales probed, within the 
limit of sensitivity discussed in Section M. Convincing ev- 
idence of thermal equilibration can be provided by event- 
by-event observables. Our data is consistent with local 
thermal equilibrium, understood as an absence of phys- 
ically distinguished scales between the scale of a hadron 
and the scale of the system, or scale invariance of fluctu- 
ations 1 33 1 ("white noise"). However to probe equilibra- 
tion directly with this method, texture sensitivity at least 
down to the typical fireball (cluster) sizes observed in 
pN collisions in cosmic rays and accelerator experiments 
p4| , |35| ] would be necessary. In the absence of such direct 
evidence, the non-observation of critical fluctuations can 
imply either absence of the second order phase transition 
or presence of thermal equilibration - the latter voids the 
criticality signature, according to Rajagopal and Wilczek 

13. 



VII. CONCLUSION 

We have developed a method of measuring the dynamic 
component of local fluctuations in charged particle den- 
sity in pseudorapidity and azimuthal angle, and applied 
the analysis to Pb + Pb collisions measured by the NA44 
experiment. Comparison of the data to a simple Monte 
Carlo texture event generator indicates that sensitivity to 
pseudorapidity density clusters of > 3% is accomplished 
in this experiment. The probability of encountering a 
real event whose dynamic azimuthal texture exceeds in 
strength that of a random mixed event by 3 RMS, is be- 
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low 10%. For the pseudorapidity texture, the respective 
probability is below 5%. 

We conclude that this novel method of event-by-event 
analysis, sensitive to particular signatures of first and 
second order phase transitions, does not reveal such sig- 
natures in 158 GeV/A Pb + Pb collisions at the SPS. 
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TABLE I: Sources of background texture (dynamic and static) 
and their treatment. The irreducible remainder estimate is 
quoted for diagonal texture correlation in the 326 < dN/ dr/ < 
398 bin, and is expressed in the units of a 2 /(dEMip) 2 ', see text 
for information on how it was obtained. 



